Automated quantification of nonperfusion in the retina using optical coherence tomography angiography

ABSTRACT

Disclosed are methods and systems for measuring areas of nonperfusion in the retina using OCT imaging. The disclosed methods and systems allow for the automated segmentation and quantification of avascular areas of the retina utilizing information obtained from both structural OCT and OCT angiography (OCTA) data. The disclosed methods include filtering approaches which enhance vessel structure while suppressing noise, dynamic thresholding approaches to mitigate the detrimental effects of within-scan variability and low scan quality, and distance transform-based approaches to improve detection of ischemic regions. When combined with methods such as projection-resolved OCTA, the sensitivity to detect nonperfusion within different plexuses of the inner retina is demonstrated. In the clinical setting of diabetic retinopathy, the disclosed methods and systems show high sensitivity and specificity to detect the mild non-proliferative form of the disease with high reproducibility.

CROSS REFERENCE TO RELATED APPLICATION

The present application claims priority to U.S. Provisional ApplicationNo. 62/364,788, which was filed on Jul. 20, 2016, and titled “AUTOMATEDQUANTIFICATION OF NONPERFUSION IN THE RETINA USING OPTICAL COHERENCETOMOGRAPHY ANGIOGRAPHY,” and which is hereby incorporated by referenceherein.

ACKNOWLEDGEMENT OF GOVERNMENT SUPPORT

Aspects of the present disclosure were made with the support of theUnited States government under the terms of grant number R01 EY023285,DP3 DK104397, and R01 EY024544 awarded by the National Institutes ofHealth. The United States government has certain rights to thisinvention.

FIELD

Generally, the field involves methods of imaging using optical coherencetomography. In particular, the field involves methods of visualizingareas of nonperfusion in the retina using optical coherence tomographyangiography.

BACKGROUND

Diabetic retinopathy (DR) is characterized by capillary nonperfusion,vascular hyperpermeability, and neovascularization, and is a leadingcause of blindness. Capillary nonperfusion, in particular, is an earlyindicator of DR and increases with severity of DR. The Early Treatmentof Diabetic Retinopathy Study (ETDRS) qualitatively evaluated macularischemia using fluorescein angiography (FA) and found it to havepredictive value for progression of disease. Such rigorous grading of FAas used in the ETDRS study, however, is impractical for clinicalpractice. FA also suffers from several additional shortcomings,including dependence on early transit for macular capillary details, dyeleakage which obscures of details of vascular structure, and variabilityof contrast. These limitations have hampered the clinicians' ability toassess nonperfusion objectively using FA. Thus, there remains an unmetneed for an accurate, objective, and automated method to evaluatemacular ischemia. Such a method would provide a valuable biomarker forDR with in-clinic applicability.

Optical coherence tomography angiography (OCTA) is an extension ofstructural optical coherence tomography (OCT) that can providehigh-contrast imaging of capillary details without the need for dyeinjection. As such, OCTA provides a method to objectively evaluatecapillary structure and health in a clinical setting. OCTA data can beused to automatically quantify total avascular area (AA) in the innerretina, and can detect DR with high sensitivity and specificity forpatients having proliferative DR, the more advanced stage of thedisease. However, detection of less severe forms of DR nonproliferativediabetic retinopathy (NPDR) remains a challenge. Part of the challengeis due to retinal anatomy: the inner retina is composed of threedistinct, but tightly stacked, plexuses that are difficult to resolveinto individual layers via OCTA because of projection artifacts.However, a new technique called Projection-resolved (PR) OCTA can beused to reduce the problem of projection artifacts blurring the plexusestogether (e.g., as described in Zhang M, Hwang T S, Campbell J P, et al.Projection-resolved optical coherence tomographic angiography.Biomedical Optics Express 2016; 7:816-828, hereby incorporated byreference herein). Using PR-OCTA, individual plexuses revealed moreareas of capillary drop out than seen in inner retinal angiograms of theplexuses. However, this increased sensitivity to detect capillarynonperfusion within individual plexuses is still susceptible tosubstantial noise artifacts. These noise artifacts degrade the integrityof imaged capillary structures and introduce spurious signal noise innon-capillary space, making quantification of AA less reliable and lessreproducible. Thus, there remains a need for improved image processingmethods to calculate AA in a robust and reliable manner.

BRIEF DESCRIPTION OF THE SEVERAL VIEWS OF THE DRAWINGS

FIG. 1 is a flowchart showing an example method for quantification ofcapillary nonperfusion, in accordance with various embodiments.

FIGS. 2A-2F are a panel of six images showing the effect of using fixedversus reflectance-adjusted thresholding on a normal eye having a lowsignal strength region (marked by circle) caused by a vitreous opacity.The low signal region seen on the original angiogram in FIG. 2A persistswith vesselness filter (shown in FIG. 2B) and binary filter with fixedthreshold (shown in FIG. 2C), falsely simulating capillary drop out.FIG. 2D shows an en face map of ganglion cell layer (GCL) and innerplexiform layer (IPL) reflectance amplitude map, which is filtered andscaled to create the reflectance-adjusted threshold image shown in FIG.2E. Using the reflectance-adjusted threshold image, the binary vesselmask (shown in FIG. 2F) does not show an area of false capillary dropout in the low reflectance area.

FIGS. 3A-3H are a panel of eight images showing an example ofnonperfusion detection in the superficial plexus of a NPDR eye. Theoriginal angiogram (shown in FIG. 3A) is enhanced by a vesselness filter(shown in FIG. 3B) and processed with binary vessel mask withreflectance-adjusted thresholding (shown in FIG. 3C), which is used tocreate a vessel distance map (shown in FIG. 3D) by applying a distancetransform. Morphologic operations are then used to remove regions withvessel distance less than 4 pixels (40 μm) (as shown in FIG. 3E), theneroded by a 5-pixel-wide square kernel and areas with minimum areasmaller than 8 pixels or minor axis smaller than 2 pixels discarded (asshown in FIG. 3F). The remaining regions were dilated by 7-pixel-widesquare kernel pixels (as shown in FIG. 3G). FIG. 3H shows the resultingavascular area (light blue) overlaid on the enhanced angiogram. It willbe apparent that other embodiments may use other pixel/distancethresholds for the operations described above with respect to FIGS.3A-3H.

FIGS. 4A-4H are a panel of eight images depicting the avascular area(AA) in two scans of a normal eye with high signal strength index (SSI),a measure of scan quality. The AA is limited to the foveal avascularzone and is similar in size and shape in all scans whether in individualplexuses or combined inner retinal angiogram. The extrafoveal AA (EAA)is defined as AA outside the white 1 mm circle (e.g., as shown in FIG.4A).

FIGS. 5A-5H are a panel of eight images depicting the AA scans of twonormal control eyes with relatively low signal strength index (SSI), ameasure of scan quality. Scans in each of the three individual plexusesand the combined inner retinal angiogram are shown. No avascular areaoutside the foveal avascular zone is detected.

FIGS. 6A-6H are a panel of eight images depicting the AA in two eyeswith mild NPDR. The three individual plexuses show incongruent areas ofcapillary nonperfusion that are not detected in the combined innerretinal angiogram.

FIGS. 7A, 7B1, 7B2, 7C1, 7C2, and 7D are a set of images showing anexample of application of a two-level distance transform applied to abinary vessel mask. The original binary vessel mask (shown in FIG. 7A)was denoised at two different denoising levels (shown in FIGS. 7B1 and7B2, respectively), then each distance transformed separately (shown inFIGS. 7C1 and 7C2, respectively), and the results averaged to generate afinal vessel distance map (shown in FIG. 7D).

FIGS. 8A1, 8A2, 8B1, 8B2, 8C1, and 8C2 are a tabular presentation ofexample image processing operations that illustrate the advantage ofusing a two-level distance transform to detect nonperfusion area.Application of a level-1 distance transform correctly detectsnonperfusion in ischemia (shown in FIG. 8A1) but falsely detectsnonperfusion in a normal eye (shown in FIG. 8A2). Application of alevel-2 distance transform correctly detects no areas of nonperfusion inthe normal eye (shown in FIG. 8B2) but fails to detect nonperfusion inischemia (shown in FIG. 8B1). Application of a two-level distancetransform with averaging correctly characterizes both the ischemia (FIG.8C1) and normal (FIG. 8C2) cases.

FIGS. 9A-9E are a set of images showing an example of thresholding andmorphological operations applied to a distance map to produce a finaldetected nonperfusion region. FIG. 9A is a vessel distance map. FIG. 9Bdepcts an initial threshold of distance map. As shown in FIG. 9C,regions with vessel distance less than 4 pixels (40 μm) are eroded by a5-pixel-wide square kernel. As shown in FIG. 9D, regions with minimumarea smaller than 8 pixels or minor axis smaller than 2 pixels arediscarded (i.e., filtered). As shown in FIG. 9E, remaining regions aredilated by 7-pixel-wide square kernel pixels, thereby generating thefinal detected nonperfusion region.

FIG. 10 schematically shows an example system for processing OCT and OCTangiography datasets to detect areas of nonperfusion in accordance withthe disclosure.

FIG. 11 schematically shows an example of a computing system inaccordance with the disclosure.

DETAILED DESCRIPTION

Disclosed are methods and systems for identifying and measuring areas ofnonperfusion in the retina. The disclosed methods and systems areapplied to OCT angiograms derived from OCT angiography datasets. Inembodiments, these OCT angiography datasets may be processed usingprojection-resolved OCTA (PR-OCTA) techniques to reduce shadowgraphicflow projection artifacts to improve resolution and quantification ofvascular and avascular areas within different retinal layers.

The disclosed methods include filtering approaches that enhance thevascular structures present in an angiogram. This enhancement involvesapplication of filter that maintains the vessel-like structure presentin the angiogram image while also denoising (removing noise) the imageto improve quantification. In embodiments, a “vesselness” filter basedon eigenvalue analysis of the local microstructure within the angiogramis applied. In further embodiments, the vesselness filter can be appliedat multiple scales to improve performance.

Also disclosed is a dynamic thresholding method that improvesclassification of flow pixels within an OCTA angiogram compared to fixed(i.e., constant) threshold approaches. In an embodiment disclosedherein, a reflectance-adjusted thresholding approach based on structuralOCT images is presented that allows for dynamic or adaptive thresholdingof angiogram images. In a disclosed embodiment, this approach producesan improved binary vessel mask segmentation for use in subsequent imageprocessing operations. An aspect of reflectance-adjusted thresholdingmethod is that it reduces sensitivity to overall scan quality, therebyimproving quantification of vascular structure and nonperfusion areaacross multiple visits and scanning conditions. This approach alsoreduces sensitivity to within-visit variation due to locally reducedsignals, which can be caused by, for example, vitreous opacity.

In an embodiment, a reflectance-adjusted thresholding image for use indynamic thresholding can be produced by generating an en facereflectance image from structural OCT data, log-transforming the en facereflectance values, and applying a smoothing filter (for example, aGaussian filter). The resultant image values may be further scaled andoffset to generate a functional form of the reflectance-adjustedthresholding relationship. The scaling and offset constants of thefunctional form can be derived, in some embodiments, from a training setof angiogram images, wherein the training angiograms are generated usingreflectance and flow data within specific layers of the retina. Amaximum threshold may also be incorporated into the reflectance-adjustedthresholding image at the location of the foveal avascular zone (FAZ) toimprove performance.

A further aspect of the disclosed methods is the use of a distancetransform to calculate a vessel distance map from a binary vessel maskfor use in quantifying avascular area. In an embodiment disclosedherein, a Euclidean distance transform is utilized. In a specificembodiment, a method for calculating a two-level vessel distance map isdisclosed, wherein a binary vessel mask is denoised at two differentlevels and a distance transform applied to each to produce two vesseldistance maps. These two vessel distance maps are then combined (forexample, by averaging) to produce a single vessel distance map for usein characterizing nonperfusion area. An aspect of such a two-levelvessel distance map is that it improves the sensitivity to detectischemia in diseased eyes while maintaining a high specificity to avoidischemia detection in healthy eyes.

Further disclosed are methods to modify the vessel distance map toimprove its performance in mapping areas of nonperfusion. In anembodiment, the vessel distance map may be thresholded into a binaryimage and then subjected to morphological operations to remove spuriousregions that arise from the threshold operation. As disclosed herein,such morphological operations may include erosion to pare the binaryimage, filtering to remove regions below a specified area or minor axislength, and/or dilation. The resultant improved vessel distance map canbe used as a nonperfusion map to highlight the locations of avascularregions in the angiogram and calculate avascular area depicted in theoriginal angiogram.

In a clinical study using the methods presented herein, the disclosedmethods are shown to have high sensitivity (able to detect mild NPDR at94.6% sensitivity with 96% specificity) and high reproducibility. It isalso demonstrated that the methods are insensitive to scan quality.

Also disclosed herein is a system for measuring nonperfusion area usingthe disclosed methods. The system may include an OCT device configuredto acquire OCT angiography data in functional connection with acomputing device having a logic subsystem and data holding capabilities.In embodiments, the computing device receives data from the OCT deviceand performs one or more of the methods described herein.

In the present detailed description, reference is made to theaccompanying drawings which form a part hereof, and in which are shownby way of illustration embodiments that can be practiced. It is to beunderstood that other embodiments can be utilized and structural orlogical changes can be made without departing from the scope. Therefore,the following detailed description is not to be taken in a limitingsense, and the scope of embodiments is defined by the appended claimsand their equivalents.

Various operations can be described as multiple discrete operations inturn, in a manner that can be helpful in understanding embodiments;however, the order of description should not be construed to imply thatthese operations are order dependent.

The description may use the terms “embodiment” or “embodiments,” whichmay each refer to one or more of the same or different embodiments.Furthermore, the terms “comprising,” “including,” “having,” and thelike, as used with respect to embodiments, are synonymous.

In various embodiments, structure and/or flow information of a samplecan be obtained using OCT (structure) and OCT angiography (flow) imagingbased on the detection of spectral interference. Such imaging can betwo-dimensional (2-D) or three-dimensional (3-D), depending on theapplication. Structural imaging can be of an extended depth rangerelative to prior methods, and flow imaging can be performed in realtime. One or both of structural imaging and flow imaging as disclosedherein can be enlisted for producing 2-D or 3-D images.

Unless otherwise noted or explained, all technical and scientific termsused herein are used according to conventional usage and have the samemeaning as commonly understood by one of ordinary skill in the art whichthe disclosure belongs. Although methods, systems, andapparatuses/materials similar or equivalent to those described hereincan be used in the practice or testing of the present disclosure,suitable methods, systems, and apparatuses/materials are describedbelow.

All publications, patent applications, patents, and other referencesmentioned herein are incorporated by reference in their entirety. Incase of conflict, the present specification, including explanation ofterms, will control. In addition, the methods, systems, apparatuses,materials, and examples are illustrative only and not intended to belimiting.

In order to facilitate review of the various embodiments of thedisclosure, the following explanation of specific terms is provided:

A-scan: A reflectivity profile that contains information about spatialdimensions and location of structures within an item of interest. AnA-scan is an axial scan directed along the optical axis of the OCTdevice and penetrates the sample being imaged. The A-scan encodesreflectivity information (for example, signal intensity) as a functionof depth (z-direction).

B-scan: A cross-sectional tomograph that can be achieved by laterallycombining a series of axial depth scans (i.e., A-scans) in thex-direction or y-direction. A B-scan encodes planar cross-sectionalinformation from the sample and is typically presented as an image.Thus, a B-scan can be called a cross sectional image.

Dataset: As used herein, a dataset is an ordered-array representation ofstored data values that encodes relative spatial location inrow-column-depth (x-y-z axes) format. In the context of OCT, as usedherein, a dataset can be conceptualized as a three dimensional array ofvoxels, each voxel having an associated value (for example, an intensityvalue or a decorrelation value). An A-scan corresponds to a set ofcollinear voxels along the depth (z-axis) direction of the dataset; aB-scan is made up of set of adjacent A-scans combined in the row orcolumn (x- or y-axis) directions. Such a B-scan can also be referred toas an image, and its constituent voxels referred to as pixels. Acollection of adjacent B-scans can be combined form a 3D volumetric setof voxel data referred to as a 3D image. In the system and methodsdescribed herein, the dataset obtained by an OCT scanning device istermed a “structural OCT” dataset whose values can, for example, becomplex numbers carrying intensity and phase information. Thisstructural OCT dataset can be used to calculate a corresponding datasettermed an “OCT angiography” dataset of decorrelation values reflectingflow within the imaged sample. There is a correspondence between thevoxels of the structural OCT dataset and the OCT angiography dataset.Thus, values from the datasets can be “overlaid” to present compositeimages of structure and flow (e.g., tissue microstructure and bloodflow) or otherwise combined or compared.

En Face angiogram: OCT angiography data can be presented as a 2Dprojection of the three-dimensional dataset onto a single planar imagecalled an en face angiogram. Construction of such an en face angiogramrequires the specification of the upper and lower depth extents thatenclose the region of interest within the retina OCT scan to beprojected onto the angiogram image. These upper and lower depth extentscan be specified as the boundaries between different layers of theretina (e.g., the voxels between the inner limiting membrane and outerplexiform layer could be used to generate an en face angiogram of theinner retina). Once generated, the en face angiogram image may be usedto quantify various features of the retinal vasculature as describedherein. This quantification typically involves the setting of athreshold value to differentiate, for example, the pixels that representactive vasculature from static tissue within the angiogram. These enface angiograms can be interpreted in a manner similar to traditionalangiography techniques such as fluorescein angiography (FA) orindocyanine green (ICG) angiography, and are thus well-suited forclinical use. It is also common to generate en face images fromstructural OCT data in a manner analogous to that used to generate enface angiograms. Angiograms from different layers may also becolor-coded and overlaid to present composite angiograms with encodeddepth information; structural en face images may also be included insuch composite image generation.

Disclosed herein is a method for automated quantification of retinalnonperfusion using OCTA in combination with several image processingtechniques. In an embodiment, the disclosed method may be combined withthe technique of projection-resolved OCT angiography (PR-OCTA) (e.g., asdescribed in U.S. patent application Ser. No. 15/374,872, herebyincorporated by reference herein), to effectively detect areas ofnonperfusion in distinct vessel networks that are normally difficult toresolve separately due to shadowgraphic projection artifacts. Forexample, as noted above, PR-OCTA allows the visualization of threedistinct retinal plexuses and detection of capillary abnormalities notvisible in inner retinal angiograms that combine these plexuses.However, even with PR-OCTA, the presence of noise artifacts within theacquired datasets introduces measurement errors that limits the reliablequantification of nonperfusion, and thus, limits the ability todistinguish patients with early or mild nonproliferative DR from normalpatients.

In an example method disclosed herein, a noisy retinal angiogramcontaining blood flow information (e.g., decorrelation values) isfiltered in such a way as to preserve the continuity capillary networkstructure in the image while also removing non-capillary noise. In anembodiment, a so-called “vesselness” filter may be used to enhance theinherent structure of the capillary network. In embodiments, a two-scalevesselness filter may be applied. Such a filter provides denoising ofthe angiogram while still enhancing both large vessels and capillaries.The denoising capability of such a filter serves to remove spuriouspixels that can interfere with the detection of actual gaps betweenvessels, an important consideration in quantification of avascular area.

The resultant enhanced angiogram may further be segmented to create abinary vessel mask for use in quantifying capillary perfusion andidentifying areas of nonperfusion. In an embodiment described herein, areflectance-adjusted thresholding approach is disclosed to create such abinary vessel mask. In various embodiments, decorrelation values of theenhanced angiogram may be compared to a threshold and given a firstbinary value if they meet the threshold and a second binary value ifthey do not meet the threshold, thereby generating the binary vesselmask. The value of the threshold may be based on the local reflectancesignal strength (e.g., for the given A-scan). The disclosedreflectance-adjusted approach reduces the sensitivity of the method tooverall scan quality and to reflectance variations within a single scan.Such variations can be introduced, for example, by vitreous opacitieswhich attenuate the optical signal and can falsely simulate capillarydropout. By using local reflectance signal strength as a reference forseparating flow signal and background noise, false positive detection ofnonperfusion is minimized and the repeatability of measurement improved.

The method may further employ a distance transform approach to processthe binary vessel mask. In embodiments, the binary vessel mask may bedenoised at multiple different levels (e.g., to remove different amountsof noise), and then a distance transform calculated for each of themultiple denoised masks to produce multiple distance maps. The resultantmultiple distance maps may be combined, for instance by a weightedaveraging approach, to produce a final distance map. An aspect of thedistance transform is that it does not give disproportionate value tolarger vessels. In addition, the use of the distance transform atmultiple denoising levels overcomes the disadvantage that the distancetransform can be sensitive to noise. Employing a multiple-level distancetransform improves the sensitivity to ischemia while maintaining a highspecificity to healthy subjects. In an embodiment described herein,application of a two-level distance map minimized disproportionatecontribution of large vessels in determining vessel density abnormalityand eliminated the undesirable detection of normal avascular area alonglarge vessels of the superficial plexus as AA. In embodiments of themethods disclosed herein, morphologic operations may be applied to thecalculated distance map to remove small regions while retaining largerregular shaped regions. Such morphological operations further reducefalse positive detection of AA.

FIG. 1 shows a flowchart of an example implementation of a method 100 inaccordance with embodiments described herein. The image processingstrategy for quantification of nonperfusion area depicted in theflowchart is composed of three principal steps: (i) pre-processing; (ii)vessel distance transform; and (iii) morphological operations. Thesesteps are described below using specific formulations of the imageprocessing filters in an example workflow. It will be understood by oneskilled in the art that other variations and formulations exist for thedifferent aspects of the method 100.

(i) Pre-Processing

In an embodiment, at 102 of the method 100, an original angiogram (e.g.,as shown in FIG. 2A and FIG. 3A) is enhanced using a vesselness filter.In some embodiments, the vesselness filter may be a two-scale vesselnessfilter, such as a two-scale (σ=1 and 2 pixels) Frangi vesselness filter.Example resultant images of an enhanced angiogram after application ofthe vesselness filter at 102 are shown in FIG. 2B and FIG. 3B. Thisfilter may enhance vessels by obtaining a vesselness measure on thebasis of eigenvalues of the second order local structure of theangiogram (Hessian) and suppresses background noise (e.g., as describedin Fraz M M, Remagnino P, Hoppe A, et al. Blood vessel segmentationmethodologies in retinal images: A survey. Computer Methods and Programsin Biomedicine 2013; 108:407-433; Sofka M, Stewart C V. Retinal VesselCenterline Extraction Using Multiscale Matched Filters, Confidence andEdge Measures. IEEE Transactions on Medical Imaging 2006; 25:1531-1546;and/or Zhang H F, Maslov K, Li M-L, Stoica G, Wang L V. In vivovolumetric imaging of subcutaneous microvasculature by photoacousticmicroscopy. Optics Express 2006; 14:9317-9323, each of which is herebyincorporated by reference herein).

At 104 of the method 100, the enhanced angiogram is thresholded togenerate a binary vessel mask. The binary vessel mask may distinguishflow signal from the background noise. In some embodiments, a fixed(i.e., constant) threshold value may be used for the thresholdingoperation at 104. However, the flow noise floor depends on the OCTreflectance signal, which may vary among scans and even within a singlescan. Consequently, with a fixed threshold cutoff, signal strengthinstability can cause within-visit variation leading to false detectionof capillary dropout from locally reduced signal (e.g., local signalattenuation caused by vitreous opacity). An example of thresholdingusing a fixed threshold value is shown in FIG. 2C.

As an alternative to fixed value thresholding, a reflectance-adjustedthreshold approach is disclosed that uses en face structural OCT of oneor more reference layers (e.g., of the ganglion cell layer (GCL) and theinner plexiform layer (IPL)). The GCL and IPL retinal layers havemoderate reflectance and are nearly free of cysts and exudates, makingthem good reference layers. In an embodiment, at 106 of the method 100,a reflectance-adjusted threshold image may be generated from areflectance image (e.g., a GCL+ICL reflectance image). For example, thereflectance-adjusted threshold image may be generated by computing themean projection of reflectance between the GCL and IPL layers (e.g., asshown in FIG. 2D), then taking the logarithm of those mean reflectionvalues, and applying a smoothing filter. In a specific embodiment, thelogarithm (S) of the mean projection is taken and filtered with a 15×15pixels Gaussian operator G, having a standard deviation of 8 pixels, andthen scaled. An example of the resultant reflectance-adjusted thresholdimage, T_(xy), is shown in FIG. 2E. Mathematically, these operations maybe represented by the following Equation (1):

T _(xy) =a×G(S)_(xy) +b  (1)

The filter size and standard deviation in Equation 1 can be setempirically. The value a is a scaling parameter, and the value b is anoffset parameter. The values a and b may be trained, for example, fromscans having in-scan signal variation. In the examples presented herein,the scaling and offset parameters a and b, respectively, were trainedusing scans from two control eyes with in-scan signal variation and twoNPDR eyes, identifying maximum capillary nonperfusion in NPDR while onlyidentifying the foveal avascular zone (FAZ) in control eyes. The valueof a was determined to be 8×10⁻⁵ and b was 8×10⁻² for the specificimages analyzed. A relatively high threshold (maximum of T) was assignedto 0.6 mm diameter central area to facilitate noise suppression in FAZ(see FIG. 2E). It should be noted that parameter values presented hereinfor the reflectance-adjusted thresholding are not intended to belimiting and one skilled in the art will understand the that they can bemodified to achieve the desired image processing effects and to improvethe quality of results.

As previously discussed, at 104 of the method 100, thereflectance-adjusted threshold image is compared with the enhancedangiogram to classify pixels as capillary or static tissue. Theresultant image is a binary vessel mask by generated by adaptive(dynamic) thresholding of the angiogram. Example resultant binary vesselmask images are shown in FIG. 2F and FIG. 3C. The adaptive thresholdingprocedure can be expressed mathematically according to Equation (2):

$\begin{matrix}{B_{xy} = \left\{ \begin{matrix}{1,} & {{{if}\mspace{14mu} D_{xy}} > T_{xy}} \\{0,} & {{{if}\mspace{14mu} D_{xy}} < T_{xy}}\end{matrix} \right.} & (2)\end{matrix}$

where B_(xy) is the binary vessel mask image, D_(xy) is the enhancedangiogram (e.g., the decorrelation values of the enhanced angiogram atthe corresponding x-y positions), and T_(xy) is the reflectance-adjustedthreshold image (e.g., the threshold values of the reflectance-adjustedthreshold image at the corresponding x-y positions).

As shown in FIG. 2, false detection of capillary nonperfusion can beavoided by using the above-described reflectance-adjusted thresholdapproach (FIG. 2F) compared to fixed thresholding (FIG. 2C).

(ii) Vessel Distance Transform

At 108 of the method 100, the binary vessel mask obtained from theprevious operations may be used to calculate a vessel distance map. Thevessel distance map may be useful for characterizing the non-vesselspace. In an embodiment, vessel distance map is calculated by applying adistance transform to the binary vessel mask image. In the examplespresented herein, a Euclidean distance transform was applied to thebinary vessel mask, and an example resultant vessel distance map isshown in FIG. 3D). The Euclidean distance transform of a binary imageassigns a number that is the distance between that pixel and the nearestnonzero pixel of the binary image (e.g., as described in [34], herebyincorporated by reference herein), and may be calculated according toEquation (3):

DT _(xy)=min_(B) _(x′y′) ₌₁√{square root over ((x−x′)²+(y−y′)²)}   (3)

where DT_(xy) is the distance transform, (x, y) is the given pixel ofthe binary vessel mask, (x′, y′) is the location of other pixels of thebinary vessel mask (e.g., to identify the nearest vessel to the givenpixel), and B_(x′y′)=1 corresponds to pixels of the binary vessel maskthat have a value of 1 (to indicate a vessel). Thus, Equation (3)minimizes the value of (x−x′)²+(y−y′)² for values of (x′, y′) for whichB_(x′y′)=1 to identify the nearest vessel to the given pixel.Accordingly, the values in the vessel distance map represent thedistance of a given pixel (x, y) to its nearest vessel (x′, y′). Anadvantage of the distance transform-based approach is that it does notgive disproportionate value to larger vessels like the vessel densitymap, allowing more accurate detection of areas of capillarynonperfusion. In further embodiments, the image to which the distancetransform is to be applied may be denoised prior to distance mapcalculation. Additionally, or alternatively, multiple distance maps canbe calculated for an image which has been denoised at different levels,and the resultant maps combined. In the examples presented herein, theEuclidean distance transform is applied to a binary vessel mask that hasbeen denoised at two different levels, and then results of the two mapsaveraged. This method ignores the normal avascular long-strip area alongthe large vessels, which is a commonly detected false positive regionwith simple thresholding methods.(iii) Morphological Operations

At 110 of the method 100, the nonperfusion map is generated from thevessel distance map. For example, the nonperfusion map may be generatedby thresholding the vessel distance map and applying one or moremorphological operations. In embodiments, the threshold level can bebased on the statistical characteristics of the computed vessel distancemap, or a value may be set empirically. In the examples presentedherein, a threshold criterion of DT>4 was used, but this level may varydepending on application, as will be understood by one skilled in theart. FIG. 3E shows an example of a vessel distance map that has beenthresholded at this DT>4 level. The morphological operations applied tothe resultant thresholded distance map can include erosion, filtering,and dilation operations, as well as other appropriate processingtechniques for processing images. For the examples presented herein,morphological operations include: (1) an erosion by 5-pixel-wide squarekernel, (2) elimination of areas smaller than 8 pixels or whose minoraxis length smaller than 2 pixels (FIG. 3F), (3) dilation by7-pixel-wide square kernel (FIG. 3G). These operations and parametervalues, however, are not intended to be limiting and one skilled in theart will understand the application of morphological operations toimprove the quality of results. The final result of these morphologicaloperations represent the final detection of nonperfusion regions in theangiogram. These results may be pseudo-colored and overlaid onto theoriginal or enhanced angiogram to aid interpretation of results, forexample, as shown in FIG. 3H.

The methods and operations described herein may be used to quantify AA,for example, in each of the three plexuses of the inner retina or in thetotal inner retinal angiogram encompassing all three plexuses. Anadvantage of the disclosed methods is that they limit detected AA to berelatively large, smooth contiguous areas and avoid detectingnon-physiologic areas as AA.

EXAMPLES Example 1—Detection of Mild NPDR Data Acquisition

Healthy volunteers and diabetic participants with mild NPDR wererecruited from Casey Eye Institute of Oregon Health and ScienceUniversity (OHSU). Eyes with non-diabetic macular pathology, mediaopacity or other significant eye disease were excluded. An informedconsent was obtained and the study was approved by the InstitutionalReview Board of OHSU. The study complied with the Declaration ofHelsinki and Health Insurance Portability and Accountability Act. Aclinical examination and masked-grading of 7-field color photographsconfirmed the severity of retinopathy.

Two 3×3 mm scans with 2 mm depth were obtained in one eye of eachparticipants within a visit using a commercial spectral domain OCTsystem (RTVue-XR; Optovue, Fremont, Calif.) with a center wavelength 840nm, a full-width half maximum bandwidth of 45 nm, and an axial scan rateof 70 kHz. In the fast transverse scanning direction, 304 axial scanswere sampled to obtain a single 3 mm B-scan. Two repeated B-scans werecaptured at a fixed position before proceeding to the next location. Atotal of 304 locations along a 3 mm distance in the slow transversedirection were sampled to form a 3D data cube. All 608 B-scans in eachdata cube were acquired in 2.9 seconds. Based on the volumetric OCTreflectance signal, the scanning software computed a signal strengthindex (SSI), which is often used as an indicator of scan quality.

Data Processing

Blood flow was detected using the split-spectrum amplitude decorrelation(SSADA) (e.g., as described in Jia Y, Tan O, TokayerJ, et al.Split-spectrum amplitude-decorrelation angiography with opticalcoherence tomography. Optics Express 2012; 20:4710-4725; Gao S S, Liu G,Huang D, Jia Y. Optimization of the split-spectrumamplitude-decorrelation angiography algorithm on a spectral opticalcoherence tomography system. Optics Letters 2015; 40:2305-2308; and/orLiu L, Jia Y, Takusagawa H L, et al. Optical coherence tomographyangiography of the peripapillary retina in glaucoma. JAMA Ophthalmology2015, each of which is hereby incorporated by reference herein).Projection artifacts were suppressed by a projection-resolved OCTA(PR-OCTA) algorithm (e.g., as described in Zhang M, Hwang T S, CampbellJ P, et al. Projection-resolved optical coherence tomographicangiography. Biomedical Optics Express 2016; 7:816-828, herebyincorporated by reference herein). OCT structural images were obtainedby averaging two repeated B-scans. The structural and angiography datawere generated simultaneously on each scan. For each scan, one x-fastscan and one y-fast scan were registered and merged through anorthogonal registration algorithm to remove motion artifacts (e.g., asdescribed in Kraus M F, Potsaid B, Mayer M A, et al. Motion correctionin optical coherence tomography volumes on a per A-scan basis usingorthogonal scan patterns. Biomedical Optics Express 2012; 3:1182-1199;and/or Kraus M F, Liu J J, Schottenhamml J, et al. Quantitative 3D-OCTmotion correction with tilt and illumination correction, robustsimilarity measure and regularization. Biomedical Optics Express 2014;5:2591-2613, both of which are hereby incorporated by reference herein).

A directional graph search algorithm identified structural boundaries(e.g., as described in Zhang M, Wang J, Pechauer A D, et al. Advancedimage processing for optical coherence tomographic angiography ofmacular diseases. Biomedical Optics Express 2015; 6:4661-4675, herebyincorporated by reference herein): the inner limiting membrane (ILM),inner plexiform layer (IPL), inner nuclear layer (INL), outer plexiformlayer (OPL), outer nuclear layer (ONL). The superficial plexus wasdefined as the slab within the inner 80% of GCC and minimum thickness of61 microns. The intermediate plexus was defined as vessels in the outer20% of GCC and inner 50% of INL. The deep plexus was defined as vesselsin the remaining slab internal to the OPL and minimum thickness of 37microns. A proportional rather than fixed segmentation scheme was chosenbecause the three plexuses merge at the edge of the foveal avascularzone (FAZ) and a fixed segmentation scheme can result in arbitrary andmeaningless differences in FAZ size between the vascular plexuses. Inthe normal eye, this scheme shows each plexus reaching the edge of theFAZ and avoids differences in the measurement of FAZ of three plexusescaused by solely by the segmentation method.

These three plexuses have been well characterized histologically innon-human primates and recently in human cadaveric eyes. An en face OCTof each plexus was generated using maximum decorrelation projection ofeach slab. A combined angiogram showing the more superficial plexusesplaced on top with different color showed the relationship between theplexuses.

The total AA (TAA) and the extrafoveal AA (EAA), defined as the AAoutside the 1 mm central circle were separately tabulated. These valueswere corrected for magnification variation associated with axial lengthvariation (e.g., as described in Huang D, Chopra V, Lu A T-H, Tan O,Francis B, Varma R. Does Optic Nerve Head Size Variation AffectCircumpapillary Retinal Nerve Fiber Layer Thickness Measurement byOptical Coherence Tomography? Disc Size and RNFL Layer ThicknessMeasurement by OCT. Investigative Ophthalmology & Visual Science 2012;53:4990-4997, hereby incorporated by reference herein).

Data Analysis

The Mann-Whitney test was used to compare the detected AA between thepatients with NPDR and the healthy controls. The diagnostic accuracy ofeach parameter was assessed by sensitivity with a fixed specificity andthe area under the receiver operating characteristic curve (AROC).Repeatability of the detected AA was assessed by evaluating the pooledstandard deviation (SD) for eyes that had two scans with a signalstrength index (SSI) greater than 54. We also compared the detected AAin different plexuses in NPDR. All statistical tests were done usingSPSS, version 20 (IBM).

The AA in healthy controls identified by the automated procedure wascompared with manual grading for detection accuracy. The manual resultswere delineated by using the free selection tool in GIMP 2.8. TheJaccard coefficient, defined as the area of intersection divided by thearea of union, measured the similarity between the automated result andmanual result. The false positive and negative errors were alsocomputed.

Results

Eyes from 13 healthy volunteers (mean [SD] age: 43[13], 10 women) and 13participants with mild NPDR (mean [SD] age: 59[8], 7 women) were imagedand one eye from each participant was randomly chosen for inclusion inthe study. The SSI of the control group ranged 59-88 and the NPDR group54-81. No eye in the study had cysts in the inner retina. In the controlgroup, the automated algorithm detected AA only inside the central 1 mmcircle in combined 1-layer or 3-layer angiograms (FIG. 4) with theexception of one eye that had an extrafoveal AA of 0.01 mm². The AA innormal controls detected by automated algorithm agreed with manualgrading with a Jaccard indices 0.85, 0.82, 0.81 for superficial,intermediate, and deep angiograms, respectively (Table 3) and wasindependent of signal strength variation between or within scans (FIG.5).

In the NPDR group, the automated algorithm identified AA that werefrequently incongruent between the three plexuses. As a result, itsometimes failed to detect any AA outside the FAZ when applied tocombined inner retinal angiogram in eyes with AA in individual plexuses(FIG. 6). AA in the superficial and deep plexuses tended to be largerthan in the intermediate plexus (Table 1).

There was no significant difference in EAA or TAA between the controland NPDR groups when the algorithm was applied to the combined 1-layerangiogram. When three plexuses were individually evaluated, the EAA andTAA in each of the three plexuses were significantly larger in the NPDRgroup compared to the control group (Table 1).

Holding the specificity at 95%, the 3-plexus showed that the sum of EAAwas the most sensitive with the best diagnostic accuracy at AROC of 0.99(Table 1). TAA had lower sensitivity and diagnostic accuracy.

Within-visit repeatability was assessed using eyes with 2 scans in thesame visit with SSI greater than 54 (Table 2). 22 eyes from controlparticipants and 17 eyes of NPDR participants had 2 scans with adequatesignal strength. The pooled SD of AA ranged from 0.000 to 0.024 mm² incontrol and from 0.035 to 0.051 mm² in NPDR (Table 2). For allparameters, the pooled SD were smaller than the population SD or thedifference between the NPDR and control group for parameter with astatistical difference between the groups (Table1).

TABLE 1 Avascular area in diabetic retinopathy and control eyes mean ±SD (mm²) Control NPDR Sensitivity, % Plexus ROI (n = 13) (n = 13) pvalue^(a) (95% CI)^(b) AROC Combined EAA 0.00 ± 0.00 0.01 ± 0.04 0.16826.9 0.62 Inner retina (8.8-56.1) TAA 0.15 ± 0.09 0.20 ± 0.10 0.121  7.70.65 (1.3-36.1) Superficial EAA 0.00 ± 0.00 0.08 ± 0.09 <0.001 91.9 0.99(65.0-98.1)  TAA 0.16 ± 0.10 0.30 ± 0.17 0.007 25.4 0.79 (6.9-55.8)Intermediate EAA 0.00 ± 0.00 0.01 ± 0.02 0.050 41.5 0.69 (17.1-70.0) TAA 0.16 ± 0.10 0.23 ± 0.10 0.040 15.4 0.71 (2.4-45.5) Deep EAA 0.00 ±0.00 0.08 ± 0.13 <0.001 78.1 0.89 (47.6-95.0)  TAA 0.15 ± 0.10 0.31 ±0.22 0.022 23.1 0.75 (5.3-53.8) Total of 3 EAA 0.00 ± 0.00 0.17 ± 0.23<0.001 94.6 0.99 plexuses (67.9-99.2)  TAA 0.46 ± 0.29 0.85 ± 0.47 0.01723.1 0.75 (5.3-53.8) ROI, region of interest; SD, standard deviation;NPDR, non-proliferative diabetic retinopathy; AROC, area under thereceiver operating curve; EAA, extrafoveal avascular area outside thecentral 1 mm circle; TAA, total avascular area in the whole scan.^(a)Using Bonferroni correction for multiple analyses, the limit offalse-positive error is 0.005. ^(b)Sensitivity to detect NPDR withspecificity held at 95% on the receiver operating curve.

TABLE 2 Within-visit repeatability of avascular area measurements PooledSD (mm²) Control NPDR Plexus ROI (n = 22^(a)) (n = 17^(b)) SuperficialEAA 0.000 0.035 TAA 0.008 0.042 Intermediate EAA 0.000 0.007 TAA 0.0070.016 Deep EAA 0.009 0.050 TAA 0.011 0.077 Total in 3 plexuses EAA 0.0040.036 TAA 0.024 0.051 ROI, region of interest; SD, standard deviation;NPDR, non-proliferative diabetic retinopathy; EAA, extrafoveal avasculararea outside the central 1 mm circle; TAA, total avascular area in thewhole scan.

TABLE 3 Agreement between automated detection and manual delineation ofthe macular avascular area in normal controls Jaccard similarity Falsepositive False negative Plexus metric error (mm²) error (mm²)Superficial 0.85 ± 0.12 0.013 ± 0.009 0.004 ± 0.003 Intermediate 0.82 ±0.16 0.017 ± 0.012 0.004 ± 0.004 Deep 0.81 ± 0.16 0.019 ± 0.019 0.004 ±0.003

Discussion

Automated quantification of nonperfusion using OCTA was improved usingthe image processing techniques described herein. PR-OCT allowedvisualization of 3 distinct retinal plexuses and detection of capillaryabnormalities that were not visible in inner retinal angiograms whichcombined the three plexuses. The vesselness filter effectively removedspurious pixels between vessels, and the use of the vessel-distance mapimproved characterization of vessel abnormality by minimizing thedisproportionate contribution of large vessels and eliminating theundesirable classification of normal avascular area along large vesselsof the superficial plexus as pathologic AA. The use of local reflectancesignal strength to adjust thresholding minimized the detection of falsepositive nonperfusion areas and improved measurement repeatability.Morphologic operations further reduced false positive detection of AA.

By applying these techniques to eyes with mild NPDR, total nonperfusionarea outside the central 1 mm circle could distinguish NPDR from controleyes with high diagnostic accuracy, even in eyes where combined innerretinal angiogram did not reveal any nonperfusion. This study confirmsthe potential of automatically quantified AA using OCTA as a biomarkerin DR, even in cases where the disease is less severe.

The FAZ is known to have significant variability in normal and diabeticpatients [36-38]. With segmentation into individual plexuses, the FAZsize can be further affected by segmentation scheme. Exclusion of thecentral 1 mm circle factored out this normal variation of FAZ andallowed good diagnostic accuracy of AA quantification.

In this cohort, the EAA was larger in the superficial and deep plexusescompared to the intermediate plexus, suggesting that nonperfusion occurslater in the intermediate plexuses in DR. A prospective study with alarger cohort can be conducted to verify this.

Example 2—Two-Level Distance Transform and Morphological Operations

FIG. 7 shows a pictorial representation of a two-level distancetransform operation to generate a vessel distance map in accordance withthe methods described herein. Starting with a binary vessel mask (FIG.7A), denoising filters are applied at two different levels to producetwo images with different amounts of noise removed (FIGS. 7B1 and 7B2).A distance transform is then applied to each of these denoised images togenerate two different vessel distance maps (FIGS. 7C1 and 7C2). Thesetwo vessel distance maps are averaged to generate a single vesseldistance map (FIG. 7D) for use in subsequent operations andcalculations.

FIG. 8 shows a pictorial representation of the improvement conferred byusing a two-level distance transformation to characterize areas ofnonperfusion in ischemic and normal conditions. As shown in row Row A,when a level 1 distance transform is applied to the binary vessel maskof an ischemic eye, the proper result is obtained. However, the samelevel 1 scheme is used in a normal eye, a false area of nonperfusion isdetected. Alternatively, as shown in Row B, when a level 2 distancetransform is applied to the binary vessel mask of ischemic eye, thetransformation fails to detect the nonperfusion area. The level 2distance transform does properly characterize the normal eye as havingno area of nonperfusion. Row C shows the effect of combing the level 1and level 2 distance transforms for these cases and averaging theresults. As shown, the combined two-level distance transform approachcorrectly characterizes both the ischemic and ischemic cases. Two-leveldistance overcomes the disadvantage that distance transform is sensitiveto noise. Thus, the two-level distance transform approach is able toimprove the sensitivity to ischemia while maintaining a high specifityto healthy subjects.

FIG. 9 shows a pictorial representation of an example set ofthresholding and morphological operations that can be used to improvethe quality of a nonperfusion map generated from a vessel distance map.Beginning with a vessel density map (FIG. 9A), the thresholdingoperation results in a binary image with a large number of spuriousregions that overestimate nonperfusion area (FIG. 9B). In this specificexample, there are numerous regions with a vessel distance less than 4pixels (40 μm). An erosion operation using a 5-pixel wide kernel resultsin a pared binary image (FIG. 9C). This pared image is further filteredto reject regions having an area less than 8 pixels or a minor axissmaller than 2 pixels (FIG. 9D). A dilation operation is applied to theremaining regions using a 7-pixel-wide square kernel, and the resultantimage (FIG. 9E) taken to be the final detected nonperfusion region.

Example 3—Optical Coherence Tomography Angiography Image ProcessingSystem

FIG. 10 schematically shows an example system 1000 for OCT imageprocessing in accordance with various embodiments. System 1000 comprisesan OCT system 1002 configured to acquire an OCT image comprising OCTinterferograms and one or more processors or computing systems 1004 thatare configured to implement the various processing routines describedherein. OCT system 1000 can comprise an OCT system suitable for OCTangiography applications, e.g., a swept source OCT system or spectraldomain OCT system.

In various embodiments, an OCT system can be adapted to allow anoperator to perform various tasks. For example, an OCT system can beadapted to allow an operator to configure and/or launch various ones ofthe herein described methods. In some embodiments, an OCT system can beadapted to generate, or cause to be generated, reports of variousinformation including, for example, reports of the results of scans runon a sample.

In embodiments of OCT systems comprising a display device, data and/orother information can be displayed for an operator. In embodiments, adisplay device can be adapted to receive an input (e.g., by a touchscreen, actuation of an icon, manipulation of an input device such as ajoystick or knob, etc.) and the input can, in some cases, becommunicated (actively and/or passively) to one or more processors. Invarious embodiments, data and/or information can be displayed, and anoperator can input information in response thereto.

In some embodiments, the above described methods and processes can betied to a computing system, including one or more computers. Inparticular, the methods and processes described herein, e.g., themethods depicted in FIGS. 1, 7, and 9 described above, can beimplemented as a computer application, computer service, computer API,computer library, and/or other computer program product.

FIG. 11 schematically shows a non-limiting computing device 1100 thatcan perform one or more of the above described methods and processes.For example, computing device 1100 can represent a processor included insystem 1000 described above, and can be operatively coupled to, incommunication with, or included in an OCT system or OCT imageacquisition apparatus. Computing device 1100 is shown in simplifiedform. It is to be understood that virtually any computer architecturecan be used without departing from the scope of this disclosure. Indifferent embodiments, computing device 1100 can take the form of amicrocomputer, an integrated computer circuit, printed circuit board(PCB), microchip, a mainframe computer, server computer, desktopcomputer, laptop computer, tablet computer, home entertainment computer,network computing device, mobile computing device, mobile communicationdevice, gaming device, etc.

Computing device 1100 includes a logic subsystem 1102 and a data-holdingsubsystem 1104. Computing device 1100 can optionally include a displaysubsystem 1106, a communication subsystem 1108, an imaging subsystem1110, and/or other components not shown in FIG. 11. Computing device1100 can also optionally include user input devices such as manuallyactuated buttons, switches, keyboards, mice, game controllers, cameras,microphones, and/or touch screens, for example.

Logic subsystem 1102 can include one or more physical devices configuredto execute one or more machine-readable instructions. For example, thelogic subsystem can be configured to execute one or more instructionsthat are part of one or more applications, services, programs, routines,libraries, objects, components, data structures, or other logicalconstructs. Such instructions can be implemented to perform a task,implement a data type, transform the state of one or more devices, orotherwise arrive at a desired result.

The logic subsystem can include one or more processors that areconfigured to execute software instructions. For example, the one ormore processors can comprise physical circuitry programmed to performvarious acts described herein. Additionally or alternatively, the logicsubsystem can include one or more hardware or firmware logic machinesconfigured to execute hardware or firmware instructions. Processors ofthe logic subsystem can be single core or multicore, and the programsexecuted thereon can be configured for parallel or distributedprocessing. The logic subsystem can optionally include individualcomponents that are distributed throughout two or more devices, whichcan be remotely located and/or configured for coordinated processing.One or more aspects of the logic subsystem can be virtualized andexecuted by remotely accessible networked computing devices configuredin a cloud computing configuration.

Data-holding subsystem 1104 can include one or more physical,non-transitory, devices configured to hold data and/or instructionsexecutable by the logic subsystem to implement the herein describedmethods and processes. When such methods and processes are implemented,the state of data-holding subsystem 1104 can be transformed (e.g., tohold different data).

Data-holding subsystem 1104 can include removable media and/or built-indevices. Data-holding subsystem 1104 can include optical memory devices(e.g., CD, DVD, HD-DVD, Blu-Ray Disc, etc.), semiconductor memorydevices (e.g., RAM, EPROM, EEPROM, etc.) and/or magnetic memory devices(e.g., hard disk drive, floppy disk drive, tape drive, MRAM, etc.),among others. Data-holding subsystem 1104 can include devices with oneor more of the following characteristics: volatile, nonvolatile,dynamic, static, read/write, read-only, random access, sequentialaccess, location addressable, file addressable, and content addressable.In some embodiments, logic subsystem 1102 and data-holding subsystem1104 can be integrated into one or more common devices, such as anapplication specific integrated circuit or a system on a chip.

FIG. 11 also shows an aspect of the data-holding subsystem in the formof removable computer-readable storage media 1112, which can be used tostore and/or transfer data and/or instructions executable to implementthe herein described methods and processes. Removable computer-readablestorage media 1112 can take the form of CDs, DVDs, HD-DVDs, Blu-RayDiscs, EEPROMs, flash memory cards, USB storage devices, and/or floppydisks, among others.

When included, display subsystem 1106 can be used to present a visualrepresentation of data held by data-holding subsystem 1104. As theherein described methods and processes change the data held by thedata-holding subsystem, and thus transform the state of the data-holdingsubsystem, the state of display subsystem 1106 can likewise betransformed to visually represent changes in the underlying data.Display subsystem 1106 can include one or more display devices utilizingvirtually any type of technology. Such display devices can be combinedwith logic subsystem 1102 and/or data-holding subsystem 1104 in a sharedenclosure, or such display devices can be peripheral display devices.

When included, communication subsystem 1108 can be configured tocommunicatively couple computing device 1100 with one or more othercomputing devices. Communication subsystem 1108 can include wired and/orwireless communication devices compatible with one or more differentcommunication protocols. As non-limiting examples, the communicationsubsystem can be configured for communication via a wireless telephonenetwork, a wireless local area network, a wired local area network, awireless wide area network, a wired wide area network, etc. In someembodiments, the communication subsystem can allow computing device 1100to send and/or receive messages to and/or from other devices via anetwork such as the Internet.

When included, imaging subsystem 1110 can be used acquire and/or processany suitable image data from various sensors or imaging devices incommunication with computing device 1100. For example, imaging subsystem1110 can be configured to acquire OCT image data, e.g., interferograms,as part of an OCT system, e.g., OCT system 1002 described above. Imagingsubsystem 1110 can be combined with logic subsystem 1102 and/ordata-holding subsystem 1104 in a shared enclosure, or such imagingsubsystems can comprise periphery imaging devices. Data received fromthe imaging subsystem can be held by data-holding subsystem 1104 and/orremovable computer-readable storage media 1112, for example.

It is to be understood that the configurations and/or approachesdescribed herein are exemplary in nature, and that these specificembodiments or examples are not to be considered in a limiting sense,because numerous variations are possible. The specific routines ormethods described herein can represent one or more of any number ofprocessing strategies. As such, various acts illustrated can beperformed in the sequence illustrated, in other sequences, in parallel,or in some cases omitted. Likewise, the order of the above-describedprocesses can be changed.

The subject matter of the present disclosure includes all novel andnonobvious combinations and subcombinations of the various processes,systems and configurations, and other features, functions, acts, and/orproperties disclosed herein, as well as any and all equivalents thereof.

1. A computer-based method of measuring nonperfusion in the retina, themethod comprising: obtaining an optical coherence tomography (OCT)angiogram; thresholding the OCT angiogram, thereby generating a binaryvessel mask; generating a vessel distance map from the binary vesselmask; and generating a nonperfusion map from the vessel distance map. 2.The method of claim 1, wherein the OCT angiogram is a projectionresolved OCT angiography (PR-OCTA) angiogram.
 3. The method of claim 1,further comprising, prior to the thresholding, applying a vesselnessfilter to the OCT angiogram to enhance the vessel structure in the OCTangiogram.
 4. The method of claim 3, wherein the vesselness filter is atwo-scale vesselness filter.
 5. The method of claim 1, whereinthresholding the OCT angiogram comprises: obtaining a structural OCT;generating an en face reflectance image from the structural OCT;filtering the en face reflectance image, thereby generating areflectance-adjusted threshold image including respectivereflectance-adjusted thresholds; and thresholding the OCT angiogramusing the respective reflectance-adjusted thresholds, thereby generatingthe binary vessel mask.
 6. The method of claim 5, wherein the structuralOCT comprises the ganglion cell layer (GCL) and inner plexiform layer(IPL).
 7. The method of claim 6, wherein the thresholding the OCTangiogram using the reflectance-adjusted thresholds is performedaccording to: $B_{xy} = \left\{ {\begin{matrix}{1,} & {{{if}\mspace{14mu} D_{xy}} > T_{xy}} \\{0,} & {{{if}\mspace{14mu} D_{xy}} < T_{xy}}\end{matrix};} \right.$ wherein B_(xy) is the binary vessel mask, D_(xy)is the OCT angiogram, and T_(xy) is the reflectance-adjusted thresholdimage.
 8. The method of claim 5, wherein the filtering the en facereflectance image comprises: calculating a logarithm of the en facereflectance image, thereby generating a log-transformed en facereflectance image; applying a smoothing filter to the log-transformed enface reflectance image, thereby generating a smoothed log-transformed enface reflectance image; scaling and offsetting the smoothedlog-transformed en face reflectance image, thereby generating a scaledreflectance-adjusted threshold image; and setting a maximum threshold ina region of the scaled reflectance-adjusted threshold image, therebygenerating the reflectance-adjusted threshold image.
 9. The method ofclaim 8, wherein the scaling and offsetting is performed according to:T _(xy) =a×G(S)_(xy) +b wherein T_(xy) is the reflectance-adjustedthreshold image, G is a Gaussian operator, S corresponds to thelog-transformed en face reflectance image, a is a scaling parameter, andb is an offset parameter.
 10. The method of claim 8, wherein thesmoothing filter is a Gaussian filter.
 11. The method of claim 10,wherein the Gaussian filter is sized 15×15 pixels with a standarddeviation of 8 pixels.
 12. The method of claim 8, wherein the region ofthe scaled reflectance-adjusted threshold image is a foveal avascularzone (FAZ) region.
 13. The method of claim 12, wherein the FAZ regioncomprises a 6 mm-diameter circular area centered on the FAZ.
 14. Themethod of claim 1, wherein generating the vessel distance map from thebinary vessel mask comprises: denoising the binary vessel mask at afirst level, thereby generating a first denoised binary mask image;denoising the binary vessel mask at a second level, thereby generating asecond denoised binary mask image; applying a distance transform to thefirst denoised binary mask image, thereby generating a level-1 vesseldistance map; applying a distance transform to the second denoisedbinary mask image, thereby generating a level-2 vessel distance map; andaveraging the level-1 vessel distance map and level-2 vessel distancemap, thereby generating the vessel distance map.
 15. The method of claim14, wherein the distance transform is a Euclidean distance transformcalculated according to:DT _(xy)=min_(B) _(x′y′) ₌₁√{square root over ((x−x′)²+(y−y′)²)} whereinDT_(xy) is the distance transform, (x, y) is a given pixel of the binaryvessel mask, (x′, y′) corresponds to other pixels of the binary vesselmask, and B_(x′y′)=1 corresponds to pixels of the binary vessel maskthat have a value to indicate a vessel.
 16. The method of claim 1,wherein generating the nonperfusion map from the vessel distance mapcomprises: thresholding the vessel distance map, thereby generating afirst thresholded image; eroding the first thresholded image, therebygenerating a second thresholded image; filtering the second thresholdedimage, thereby generating a third thresholded image; and dilating thethird thresholded image, thereby generating the nonperfusion map. 17.The method of claim 16, wherein the thresholding the vessel distance mapis performed using a value of DT>4 when the vessel distance map iscalculated according to:DT _(xy)=min_(B) _(x′y′) ₌₁√{square root over ((x−x′)²+(y−y′)²)} whereinDT_(xy) is the distance transform, (x, y) is a given pixel of the binaryvessel mask, (x′, y′) corresponds to other pixels of the binary vesselmask, and B_(x′y′)=1 corresponds to pixels of the binary vessel maskthat have a value to indicate a vessel.
 18. The method of claim 16,wherein the eroding the first thresholded image is performed using a5-pixel wide square kernel.
 19. The method of claim 16, wherein thefiltering the second thresholded image comprises eliminating areas ofthe second thresholded image smaller than 8 pixels.
 20. The method ofclaim 19, wherein the filtering the second thresholded image furthercomprises eliminating areas of the second thresholded image whose minoraxis length is smaller than 2 pixels.
 21. The method of claim 16,wherein the dilating the third thresholded image is performed using a7-pixel wide square kernel.
 22. A system for quantifying nonperfusion inthe retina, comprising: an OCT system to acquire an OCT angiogram of asample; a logic subsystem; and a data holding subsystem comprisingmachine-readable instructions stored thereon that are executable by thelogic subsystem to: threshold the OCT angiogram, thereby generating abinary vessel mask; generate a vessel distance map from the binaryvessel mask; and generate a nonperfusion map from the vessel distancemap.
 23. The system of claim 22, wherein the OCT angiogram is aprojection resolved OCT angiography (PR-OCTA) angiogram.
 24. The systemof claim 22, wherein the instructions are further executable by thelogic subsystem to, prior to performing the thresholding, apply avesselness filter to the OCT angiogram to enhance the vessel structurein the OCT angiogram.
 25. The system of claim 22, wherein, to thresholdthe OCT angiogram, the logic subsystem is to: obtain a structural OCT ofthe sample; generate an en face reflectance image from the structuralOCT; filter the en face reflectance image, thereby generating areflectance-adjusted threshold image including respectivereflectance-adjusted thresholds; and threshold the OCT angiogram usingthe respective reflectance-adjusted thresholds, thereby generating thebinary vessel mask.
 26. The system of claim 25, wherein the structuralOCT comprises the ganglion cell layer (GCL) and inner plexiform layer(IPL).
 27. The system of claim 6, wherein, the logic subsystem is tothreshold the OCT angiogram using the reflectance-adjusted thresholdsaccording to: $B_{xy} = \left\{ {\begin{matrix}{1,} & {{{if}\mspace{14mu} D_{xy}} > T_{xy}} \\{0,} & {{{if}\mspace{14mu} D_{xy}} < T_{xy}}\end{matrix};} \right.$ wherein B_(xy) is the binary vessel mask, D_(xy)is the OCT angiogram, and T_(xy) is the reflectance-adjusted thresholdimage.
 28. The system of claim 25, wherein, to filter the en facereflectance image, the logic subsystem is to: calculate a logarithm ofthe en face reflectance image, thereby generating a log-transformed enface reflectance image; apply a smoothing filter to the log-transformeden face reflectance image, thereby generating a smoothed log-transformeden face reflectance image; scale and offset the smoothed log-transformeden face reflectance image, thereby generating a scaledreflectance-adjusted threshold image; and set a maximum threshold in aregion of the scaled reflectance-adjusted threshold image, therebygenerating the reflectance-adjusted threshold image.
 29. The system ofclaim 28, wherein the logic subsystem is to scale and offset thesmoothed log-transformed en face reflectance image according to:T _(xy) =a×G(S)_(xy) +b wherein T_(xy) is the reflectance-adjustedthreshold image, G is a Gaussian operator, S corresponds to thelog-transformed en face reflectance image, a is a scaling parameter, andb is an offset parameter.
 30. The system of claim 28, wherein thesmoothing filter is a Gaussian filter sized 15×15 pixels with a standarddeviation of 8 pixels.
 31. The system of claim 28, wherein the region ofthe scaled reflectance-adjusted threshold image is a foveal avascularzone (FAZ) region.
 32. The system of claim 31, wherein the FAZ regioncomprises a 6 mm-diameter circular area centered on the FAZ.
 33. Thesystem of claim 22, wherein, to generate the vessel distance map fromthe binary vessel mask, the logic subsystem is to: denoise the binaryvessel mask at a first level, thereby generating a first denoised binarymask image; denoise the binary vessel mask at a second level, therebygenerating a second denoised binary mask image; apply a distancetransform to the first denoised binary mask image, thereby generating alevel-1 vessel distance map; apply a distance transform to the seconddenoised binary mask image, thereby generating a level-2 vessel distancemap; and average the level-1 vessel distance map and level-2 vesseldistance map, thereby generating the vessel distance map.
 34. The systemof claim 33, wherein the distance transform is a Euclidean distancetransform calculated according to:DT _(xy)=min_(B) _(x′y′) ₌₁√{square root over ((x−x′)²+(y−y′)²)} whereinDT_(xy) is the distance transform, (x, y) is a given pixel of the binaryvessel mask, (x′, y′) corresponds to other pixels of the binary vesselmask, and B_(x′y′)=1 corresponds to pixels of the binary vessel maskthat have a value to indicate a vessel.
 35. The system of claim 22,wherein, to generate the nonperfusion map from the vessel distance map,the logic subsystem is to: threshold the vessel distance map, therebygenerating a first thresholded image; erode the first thresholded image,thereby generating a second thresholded image; filter the secondthresholded image, thereby generating a third thresholded image; anddilate the third thresholded image, thereby generating the nonperfusionmap.
 36. The system of claim 35, wherein the logic subsystem is tothreshold the vessel distance map using distance transforms of thevessel distance map having a value greater than 4, and wherein thevessel distance map is calculated according to:DT _(xy)=min_(B) _(x′y′) ₌₁√{square root over ((x−x′)²+(y−y′)²)} whereinDT_(xy) is the respective distance transform of the vessel distance map,(x, y) is a given pixel of the binary vessel mask, (x′, y′) correspondsto other pixels of the binary vessel mask, and B_(x′y′)=1 corresponds topixels of the binary vessel mask that have a value to indicate a vessel.37. The system of claim 36, wherein the first thresholded image iseroded using a 5-pixel wide square kernel; wherein the secondthresholded image is filtered by eliminating areas of the secondthresholded image smaller than 8 pixels and eliminating areas of thesecond thresholded image whose minor axis length is smaller than 2pixels; and wherein the third thresholded image is dilated using a7-pixel wide square kernel.